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Motivated by the indication of a new critical theory for the spin- 1/2 Heisenberg model with a 
spatially staggered anisotropy on the square lattice as suggested in pj|, we re- investigate the phase 
transition of this model induced by dimerization using first principle Monte Carlo simulations. We 
focus on studying the finite-size scaling of p a \2L and p 3 22L, where L stands for the spatial box size 
used in the simulations and p a i with i £ {1,2} is the spin-stiffness in the i-direction. Remarkably, 
while we do observe a large correction to scaling for the observable p a \2L as proposed in @], the 
data for p a 22L exhibit a good scaling behavior without any indication of a large correction. As a 
consequence, we are able to obtain a numerical value for the critical exponent v which is consistent 
with the known 0(3) result with moderate computational effort. Specifically, the numerical value 
of v we determine by fitting the data points of p a 22L to their expected scaling form is given by 
v — 0.7120(16), which agrees quantitatively with the most accurate known Monte Carlo 0(3) result 
v — 0.7112(5). Finally, while we can also obtain a result of v from the observable second Binder 
ratio Q2 which is consistent with v = 0.7112(5), the uncertainty of v calculated from Q2 is more 
than twice as large as that of v determined from p a 22L. 
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Introduction. — Heisenberg-type models have been 
one of the central research topics in condensed matter 
physics during the last two decades. The reason that 
these models have triggered great theoretical interest is 
twofold. First of all, Heisenberg-type models are rele- 
vant to real materials. Specifically, the spin- 1/2 Heisen- 
berg model on the square lattice is the appropriate model 
for understanding the undoped precursors of high T c 
cupratcs (undoped antiferromagnets) . Second, because 
of the availability of efficient Monte Carlo algorithms 
as well as the increased power of computing resources, 
properties of undoped antiferromagnets on geometrically 
non-frustrated lattices can be simulated with unprece- 
dented accuracy Hrl- As a consequence, these mod- 
els are particular suitable for examining theoretical pre- 
dictions and exploring ideas. For instance, Heisenberg- 
type models are often used to examine field theory pre- 
dictions regarding the universality class of a second or- 
der phase transition P, 0, @ • Furthermore, a new pro- 
posal of determining the low-energy constant, namely the 
spinwave velocity c of antiferromagnets with 0(2) and 
0(3) symmetry, through the squares of temporal and 
spatial winding numbers was verified to be valid and 
this idea has greatly improved the accuracy of the re- 
lated low-energy constants [HI EJ- On the one hand, 
Heisenberg-type models on geometrically non-frustrated 
lattices arc among the best quantitatively understood 
condensed matter physics systems; on the other hand, 
despite being well studied, several recent numerical inves- 
tigation of spatially anisotropic Heisenberg models have 
led to unexpected results 0, [TH, [l]| . In particular, Monte 
Carlo evidence indicates that the anisotropic Heisenberg 
model with staggered arrangement of the antifcrromag- 
netic couplings may belong to a new universality class, 
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in contradiction to the theoretical 0(3) universality pre- 
diction [l|. For example, while the most accurate Monte 
Carlo value for the critical exponent v in the 0(3) uni- 
versality class is given by v = 0.7112(5) [l4|, the cor- 
responding v determined in [ij is v = 0.689(5). Al- 
though the subtlety of calculating the critical exponent 
v from performing finite-size scaling analysis has been 
demonstrated for a similar anisotropic Heisenberg model 
on the honeycomb lattice [l5j], the discrepancy between 
v = 0.689(5) and v = 0.7112(5) observed in [J G3] re- 
mains to be understood. 

In order to clarify this issue further, several efforts have 
been devoted to studying the phase transition of this 
model induced by dimerization. For instance, an uncon- 
ventional finite-size scaling is proposed in [ly]. Further, 
it is argued that there is a large correction to scaling 
for this phase transition which leads to the unexpected 
v = 0.689(5) obtained in pj. Still, direct numerical ev- 
idence to solve this puzzle is not available yet. In this 
study we undertake the challenge of determining the 
critical exponent v by simulating the spin-1/2 Heisenberg 
model with a spatially staggered anisotropy on the square 
lattice. The relevant observables considered in this study 
for calculating the critical exponent v arc p s \2L, p S 22L 
and Qi- Here p S i with i € {1,2} are the spin stiffness in 
the i-direction, L is the box size used in the simulations 
and Q2 is the second Binder ratio which will be defined 
later. Further, we analyze in more detail the finite-size 
scaling of p s \2L and p S 22L. The reason that p s \2,L and 
p s -{2L are chosen is twofold. First of all, these two ob- 
servables can be calculated to a very high accuracy using 
loop algorithms Q. Second, one can measure p s \ and p s i 
separately. On isotropic systems, one would naturally 
use p s , which is the average of p s \ and p S 2, for the data 
analysis. However, for the anisotropic model considered 
here, we find it is useful to analyze both the data of p s \ 
and p S 2 because such a study may reveal the impact of 
anisotropy on the system. Surprisingly, as we will show 
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FIG. 1: The spatially anisotropic Heisenberg model consid- 
ered in this study. 

later, the observable p s i2L receives a much less severe 
correction than p s \2L. Hence p S 2^L is a better suited 
quantity than p s \2L for the finite-size scaling analysis. 
Indeed, with moderate computational effort, we can ob- 
tain a numerical value for v consistent with the most 
accurate 0(3) Monte Carlo result v = 0.7112(5). 

This paper is organized as follows. First, the 
anisotropic Heisenberg model and the relevant observ- 
ables studied in this work are briefly described after 
which we present our numerical results. In particular, 
the corresponding critical point as well as the critical ex- 
ponent v are determined by fitting the numerical data to 
their predicted critical behavior near the transition. A 
final section then concludes our study. 

Microscopic Model and Corresponding 
Observables. — The Heisenberg model considered 
in this study is defined by the Hamilton operator 

H = JS X ■ S y + J 'S*' ■ 3c. (!) 

{xy) (x'y') 

where J' and J are antiferromagnetic exchange couplings 
connecting nearest neighbor spins (xy) and (x'y 1 ), re- 
spectively. Figure 1 illustrates the Heisenberg model de- 
scribed by Eq. ([1]). To study the critical behavior of this 
anisotropic Heisenberg model near the transition driven 
by the anisotropy, in particular, to determine the critical 
point as well as the critical exponent v, the spin stiff- 
nesses in the 1- and 2-directions which are defined by 

Psi = ^2<W?>. ( 2 ) 

are measured in our simulations. Here /3 is the inverse 
temperature and L again refers to the spatial box size. 
Further (Wf) with i £ {1,2} is the winding number 
squared in the i-direction. In addition, the second Binder 
ratio Q2, which is defined by 



is also measured in our simulations as well. Here m z s is 
the z-component of the staggered magnetization m s = 
T? 2a;( — l) Xl+X2 &x- By carefully investigating the spa- 
tial volume and the J 1 / J dependence of p S {L as well as 



Q2, one can determine the critical point and the critical 
exponent v with high precision. 
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TABLE I: The numerical values for v and (J'/J) c calculated 
from p a 22L, p a \2L and Q2- All results are obtained using a 
second order Taylor series expansion of Eq. @ , except those 
with a star, which are determined using an expansion to third 
order. The confluent correction to is included in the fit explic- 
itly only for p a \2L. 

Determination of the Critical Point and the 
Critical Exponent v. — To study the quantum phase 
transition, we have carried out large scale Monte Carlo 
simulations using a loop algorithm. Further, to calculate 
the relevant critical exponent v and to determine the lo- 
cation of the critical point in the parameter space J' /J, 
we have employed the technique of finite-size scaling for 
certain observables. For example, if the transition is sec- 
ond order, then near the transition the observable p S i2L 
for i £ {1,2} and Q2 should be described well by the 
following finite-size scaling ansatz 

L (t) = (1 + bL-^goitL^^jS/L*), (4) 

where Ol stands for p S i2L with i £ {1,2} or Q2, t = 
(ic — j)/jc with j = (J'/ J), b is some constant, v is the 
critical exponent corresponding to the correlation length 
£, lu is the confluent correction exponent and z is the 
dynamical critical exponent which is 1 for the transition 
considered here. Finally, go is a smooth function of the 
variables tL 1 ^ and f3/L z . From Eq. @, one concludes 
that the curves for Ol corresponding to different L, as 
functions of J' / J, should intersect at the critical point 
(J'/J) c for large L. To calculate the critical exponent 
v and the critical point (J'/J) c , in the following we will 
apply the finite-size scaling formula, Eq. @, to p s \2L, 
p S 22L as well as Q2- Without loss of generality, we have 
fixed J = 1 in our simulations and have varied J' . Addi- 
tionally, the box size used in the simulations ranges from 
L = 16 to L = 136. To reach a lattice size as large as pos- 
sible, we use (3J = 2L for each L in our simulation. As 
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a result, the temperature dependence in Eq. (|4]) drops 
out. Figure [5] shows the Monte Carlo data for p s i2L, 
p S 2%L and Q2 as functions of J' /J. The figure clearly 
indicates that the phase transition is most likely second 
order since for all the observables p s \2L, p S 22L and Q2, 
the curves of different L tend to intersect near a particu- 
lar point in the parameter space J'/ J. The most striking 
observation from our results is that the observable p s \2L 
receives a much more severe correction than p S 22L. This 
can be understood from the trend of the crossing among 
these curves for different L near the transition (figure 
[3]). Therefore one expects that a better determination 
of v can be obtained by applying the finite-size scaling 
ansatz Eq. ([J| to p S 22L. Before presenting our results, 
we would like to point out that data from large volumes 
are essential in order to determine the critical exponent 
v accurately as was emphasized in [lj|. We will use the 
strategy employed in [15( for our data analysis as well. 
Let us first focus on p S 22L since this observable shows a 
good scaling behavior. Notice from figure [31 the trend of 
crossing for different L of p S 2%L indicates that the con- 
fluent correction is negligible for lattices of larger size. 
Therefore one expects that a result consistent with the 
theoretical prediction can be reached with b = in for- 
mula Q if data from large L are taken into account in 
the fit. Indeed with a Taylor expansion of Eq. Q up to 
second order in tL l ' v as well as letting b = in Eq. 
we arrive at v = 0.7120(16) and (J'/J) c = 2.51950(3) 
using the data of p S 2%L with 48 < L < 136. In obtain- 
ing the results v = 0.7120(16) and (J'/J) c = 2.51950(3), 
we have performed bootstrap sampling on the raw data 
and have carried out a large number (around 1000) of 
fits. The inclusion of higher order terms in the Taylor 
series and eliminating data of smaller L in the fits leads 
to compatible (and consistent) results with what we have 
just obtained. Notice that the v we obtain is consistent 
with the most accurate Monte Carlo result v = 0.7112(5) 
in the 0(3) universality class. Further, the critical point 
(J' I J) c = 2.51950(3) we calculate agrees with the known 
results in the literature [l|, [l6| as well. 

Having determined v = 0.7120(16) using the observ- 
able p S 22L, we turn to the calculations of v based on the 
observables p 8 \2L and Q2- First of all, we would like to 
reproduce the unexpected result v = 0.689(5) found in 
[J. Indeed, using the Monte Carlo data of p s \2L with 
L ranging from L = 24 to L = 80 (the size of L = 80 
is similar to the largest lattice (L = 72) used in [l| in 
obtaining v = 0.689(5)), we arrive at v = 0.689(3) and 
(J' I J) c = 2.5194(4), both of which are statistically con- 
sistent with those determined in [l|. Further, a numeri- 
cal value for v consistent with its theoretical prediction 
v = 0.7112(5) could never have been obtained using the 
available data for p s \2L . This implies that the correc- 
tion to scaling for p s \2L is large and the lattice data for 
p s \2L with larger L are required to reach a numerical 
value of v consistent with the theoretical expectation. 
Finally, for the observable Q2, using the leading finite- 
size scaling ansatz (i.e. letting b — in Eq. we are 
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FIG. 2: Monte Carlo data of p sl 2L (top), p a2 2L (middle), 
and Q2 (bottom). 



able to reach a value for v which agrees even quantita- 
tively with v = 0.7112(5). However, the uncertainty of v 
calculated from Q2 is more than twice as large as that of 
v determined from p S 22L. 

The results for v and (J'/J) c calculated from our 
finite-size scaling analysis are summarized in table 1. 

Conclusions. — In this paper, we re-investigated the 
critical behavior at the phase transition induced by 
dimerization of the spin- 1/2 Hciscnbcrg model with a 
spatially staggered anisotropy. Unlike the scenario sug- 
gested in [l| that an unconventional universality class 
is observed, we conclude that indeed this second order 
phase transition is well described by the 0(3) universality 
class. Our observation of p S 22L being a good observable 
for determining the critical exponent v is crucial for read- 
ing this conclusion by confirming the 0(3) critical expo- 
nent for this phase transition with high precision. While 
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FIG. 3: Crossing of p a \1L (top), p S 22L (middle), and Qi 
(bottom) for different L near the transition. 



we do observed a large correction to scaling for the ob- 
servable p s \2L as proposed in @|, the data points of p S 2^L 
show good scaling behavior. Specifically, with p S 22L. we 
can easily reach a highly accurate numerical value for v 
consistent with the theoretical predictions without tak- 
ing the confluent correction into account in the fit. The 
large correction to scaling observed for p s \lL, in principle 
should influence all obscrvablcs. Hence the most reason- 
able explanation for the good scaling behavior of p S 2^L 
shown here is that the prefactor b in Eq. Q for p S 2^L 
is very small. As a result, we are able to determine the 
expected numerical value for v using data of p S 2^L with 
moderate lattice sizes. Still, a more rigorous theoretical 
study such as investigating whether there exists a symme- 
try that protects p S 2^L from being affected by the large 
correction to scaling as suggested in [2[ will be an inter- 
esting topic to explore. For example, in @] it is argued 
that the enhanced correction to scaling observed for this 
phase transition might be due to a cubic irrelevant term 
which contains one-derivative in the 1-direction. The first 
thing one would like to understand is whether the feature 
of this irrelevant term, namely it contains one-derivative 
in the 1-direction, will lead to our observation that the 
large correction to scaling has little impact on p S 22L. 
In summary, here we present convincing numerical evi- 
dence to support that the phase transition considered in 
this study is well described by the 0(3) universality class 
prediction, at least for the critical exponent v which is 
investigated in detail in this study. Finally, whether the 
good scaling of p S 2%L observed here is a coincidence or 
is generally applicable for quantum Hciscnbcrg models 
with a similar spatially anisotropic pattern, remains an 
interesting topic for further investigation. 
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